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Abstract 

An aeroelastic analysis system for flutter and forced response analysis of turbomachines based on 
a two-dimensional linearized unsteady Euler solver has been developed. The ASTROP2 code, 
an aeroelastic stability analysis program for turbomachinery, was used as a basis for this 
development. The ASTROP2 code uses strip theory to couple a two dimensional aerodynamic 
model with a three dimensional structural model. The code was modified to include forced 
response capability. The formulation was also modified to include aeroelastic analysis with 
mistuning. A linearized unsteady Euler solver, LINFLX2D is added to model the unsteady 
aerodynamics in ASTROP2. By calculating the unsteady aerodynamic loads using LINFLX2D, 
it is possible to include the effects of transonic flow on flutter and forced response in the 
analysis. The stability is inferred from an eigenvalue analysis. The revised code, ASTROP2-LE 
for ASTROP2 code using Linearized Euler aerodynamics, is validated by comparing the 
predictions with those obtained using linear unsteady aerodynamic solutions. 


Introduction 

The aeroelastic research program at NASA Glenn Research Center is focused on flutter, and 
forced response analysis of turbomachinery. An overview of this research was presented in 
Ref. 1 . The review showed that a range of aerodynamic and structural models have been used to 
obtain the aeroelastic equations. Both time and frequency domain methods have been used to 
obtain unsteady aerodynamic forces and to solve the aeroelastic equations. It was noted that time 
domain methods require large computational time compared to frequency domain methods, and 
should only be used when non-linearities are expected, and when the need justifies the cost. 

Two approaches were used in obtaining the unsteady aerodynamic forces using frequency 
domain methods. In the first approach. Ref. 2, the unsteady aerodynamic equations are 
linearized about a uniform steady flow, there by neglecting the effects of airfoil shape, angle o 
attack and thickness. The unsteady aerodynamic models developed in Refs. 3-6 based on this 
approach are used in Refs. 7-8 to study the flutter and forced response analysis of a compressor 
rotor using a typical section structural model. Some of these models were also integrated with a 
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three dimensional structural model using strip theory in Ref. 9. However, methods developed by 
this approach are restricted to shock-free flows and lightly loaded blade rows. 


In the second approach, Ref. 10, the unsteady flow is regarded as a small amplitude perturbation 
about a non-uniform steady flow. The unsteady non-linear aerodynamic equations are linearized 
about the non-uniform steady flow, resulting in variable coefficient linear unsteady aerodynamic 
equations, which include the effects of steady aerodynamic loading due to airfoil shape, 
thickness and angle of attack. Following the second approach, Refs. 11-13 developed a nonlinear 
steady and linear unsteady aerodynamic model based on the potential equation. This unsteady 
aerodynamic model was used to study the effect of steady aerodynamic loading on flutter 
stability using a typical section structural model in Ref. 14. Subsequently this aerodynamic 
model was integrated with a three dimensional structural model using strip theory in Ref. 15. 

However, the formulation based on the potential equation requires corrections for entropy and 
flow rotation. The Euler equations can be used to correctly model rotational and entropy effects 
associated with strong shocks. Unsteady linearized Euler aerodynamic models that include the 
effect of steady aerodynamic loading were developed in Refs. 16-18. Recently, a two 
dimensional linearized Euler code, LINFLX2D, was developed under a NASA contract’ and has 
been reported in reference 19. This code is based on the non-linear Euler solver developed in 
Ref. 20. The LINFLX2D code can be used with any aeroelastic code that uses a typical section 
aeroelastic model as its basis. In Ref. 21, the unsteady aerodynamic calculations from 
LINFLX2D were coupled with MISER (Ref. 7); a frequency domain aeroelastic stability and 
response code based on a typical section structural model. Flutter and forced response 
calculations were presented for cascades in subsonic and transonic flow, with and without 
mistuning. 


The aeroelastic formulation in MISER does not represent the behavior of a three dimensional 
structure. An ideal analysis will be to couple a three dimensional structural analysis with a three 
dimensional aerodynamic analysis. This may be computationally intensive. An intermediate 
approach that couples a two dimensional aerodynamic analysis with a three dimensional 
structural model using strip theory, quasi-3D approach, may be less expensive computationally. 
This approach will provide accurate results except where three dimensional effects dominate. In 
turbomachines, where the compressors and turbines are enclosed, three dimensional effects may 
not be strong at least away from blade tip, and on stators, the strip theory can be used. The 

ASTROP2 code reported in Ref. 9 uses the strip theory approach for aeroelastic stability 
analysis. 


The primary objective of the present study is to develop a quasi-3D aeroelastic code by coupling 
2-D hnearized Euler analysis code, LINFLX2D, with a three dimensional structural model 
using strip theory. This effort uses the ASTROP2 code of Ref. 9 as the basis. The ASTROP2 
code is a frequency domain stability analysis program, and as mentioned before uses strip theory 
to couple a two dimensional aerodynamic model with a three dimensional structural model. The 
present version of ASTROP2 has provision for flutter calculations, but not for forced response. 
During the course of the research effort, the ASTROP2 code is extended to include forced 
response calculation. Even though the formulation in ASTROP2 allowed for frequency mistuned 
cascade analysis, to the authors’ knowledge, the code has never been validated and used for 
mistuned cascade analysis. In the present effort, the original formulation is modified to 
include both frequency mistuning and mode shape mistuning. The resulting code is named 
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ASTROP2-LE for ASTROP2 code using Linearized Euler aerodynamics, and validated by 
applying it to a representative blade. Brief descriptions of the formulation and method of 
analysis are given in the next section, followed by results and discussion. 


Formulation 

The aeroelastic formulation using the linearized approach requires solutions from two 
aerodynamic codes and a solution from a structural dynamic analysis code. The steady 
aerodynamic loads are obtained from a non-linear Euler solver, NPHASE, and the unsteady 
aerodynamic loads are obtained from LINFLX2D. The structural dynamic solution can be 
obtained from any finite element code, analytical solutions or measured data. The salient 
features of both the aerodynamic codes, and the aeroelastic formulation are described below. 

Aerodynamic model 

Non-linear Steady Euler Solver, NP HASE 

The steady aerodynamic model is based on the unsteady, two-dimensional, Euler equations. The 
equations in conservative differential form are solved in a time-dependent body-fitted curvilinear 
reference frame. This transformation process and the ensuing numerical method are presented in 
detail in Ref. 20. The equations are discretized and solved using a finite volume method with a 
combination of flux difference splitting and flux vector splitting scheme. In addition, limiters are 
used to control dispersive errors commonly encountered with higher order schemes. The steady 
solutions presented herein are obtained using the implicit scheme developed in Ref. 20, which is 
third order accurate in space and second order accurate in time. 

Linear Unsteady Euler Solver. LINFLX2D 

To obtain the linearized unsteady Euler equations, the dependent variables in the unsteady non- 
linear Euler equations are expanded in an asymptotic series of the form 

U = U(x) + u(x(x,t),t) + higher order terms ( 1 ) 

where, the term U(x) is of order one and the second term is of the order e. Substituting the 
expansion of Eq. 1 in the nonlinear unsteady Euler equations, and equating terms of like power in 
e, and neglecting terms of second order in £, nonlinear steady equations and linear variable 
coefficient unsteady equations are obtained. The unsteady linear equations are further simplified 
by assuming the unsteady excitations and responses are harmonic in time as 

U = U(x) + Re[u(x)exp(i cut)] (2) 

For harmonic blade motions with constant phase angle between adjacent blades (interblade phase 
angle), the values of interblade phase angle (a) that can occur are given as in Ref. 22. 

a r = 27tr/N. r = 0,1,2, ,N-1 (2) 
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N is the number of blades in the cascade. In a time domain approach, the number of blocks 
required depends on the interblade phase angle, and small phase angles may require large number 
of blocks to calculate the unsteady aerodynamic forces. However, with the linear approach, the 
periodic conditions are applied on a single extended blade passage region i.e., a region of angular 
pitch, 


6=2n/N (4) 

In solving the linear unsteady equations, the independent variables are regarded as pseudo time 
dependent. This allows solutions to be determined using conventional time -marching 
algorithms to converge the steady and the complex amplitudes of the unsteady conservation 
variables to their steady state values. For more details, see reference 19. 

Aeroelastic Model 

As mentioned before, the ASTROP2 code of Ref. 9 is selected to couple the unsteady linearized 
aerodynamic solver, LINFLX2D. The ASTROP2 code uses strip theory to couple a two 
dimensional aerodynamic model with a three dimensional structural model. The present version 
of ASTROP2 can solve only for flutter and not for forced response. For clarity and 
completeness, the aeroelastic formulation in ASTROP2 is given here with extension to include 
forced response calculation. 

ASTROP2 uses the normal mode approach for aeroelastic analysis. The equations of motion for 
the k blade of the cascade for ‘P’ normal modes can be written as 

m k Px p{m k Px , +[K? RxP {n} k Pxl = {/*&, + {f a } k Px , ( 5 ) 

where [Mf and [K] k are generalized mass and stiffness matrices, {r|} k is the generalized 
displacement vector, and is the motion dependent aerodynamic load vector, and {/ 0 }*is 

the motion independent aerodynamic load vector. The motion dependent forces cause flutter, 
and motion independent forces cause forced response (forced vibration). The elements of the 
matrices [M] and [K] are given by a free vibration analysis. The expressions for {f n °y and 
y a } using strip theory are given below. The expressions are developed in terms of the 
generalized displacement coordinates, {r| }, (Ref. 23), instead of interblade phase angle modes as 
was done in Refs. 7 and 8. 

In ASTROP2 the blades are divided into strips where the aerodynamic forces are calculated, see 
Fig. la. Each strip has two degrees of freedom, a plunging displacement, h, motion 
perpendicular to chord, and a pitching (torsion) displacement a, rotation about the leading edge 
of the strip (Fig. lb). Using the normal modal values obtained from a free vibration analysis, the 
equivalent h and a for m lh strip of the k th blade are given as summation of normal modes as 


NASA/TM — 2002-2 1 1499 


4 


( 6 . 1 ) 


K}‘ = 


1 

fc 

-s: 

XA2" 

■Kr ' 

KL 

r 

R 

pf 

Tt 

1 


J2 xP 


n, 

v 2 

n P 


I Pa l 


or 

K}‘=MW 


(6.2) 


where [(p m ] is the modal matrix for the m lh strip 


It should be noted that when using ASTROP2, only the location of the strips has to be input to 
the code. ASTROP2 calculates the h and a values at these strips, and uses in the analysis i.e. the 
user need not input the h and a values at each strip for each mode. 


Assembling for M strips, the modal matrix for the k'" blade is given by 
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The motion dependent aerodynamic forces {/''«} are given by 


2MxP 


(7) 


where 0) is the assumed frequency at which the aerodynamic forces are calculated, and{F} is 
given as 




n=l 


2MxP L J P.vl 


( 8 ) 


where N is the number of blades, and 
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Here , & m , 4, , p m are respectively the semichord, reduced frequency based on b m , length of 

the strip, and air density at the m strip. C Fq and C Fa are the lift coefficients due to unit 

amplitudes of plunging and pitching displacements, and C r Mq and C Ma are the moment 

coefficients due to unit amplitudes of plunging and pitching displacements for the r* interblade 
phase angle, respectively. 


The motion independent aerodynamic forces, {f a }* are given by 


i* 

\ Px\ 


( 11 ) 


where 
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( 12 . 2 ) 
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Here C and C r Mw are the lift and moment coefficients due to unit amplitudes of the wake, and 
w r is the amplitude of the velocity of the sinusoidal wake in the r 1 * 1 interblade phase angle mode. 
Since the forced response problem is typically for a given interblade phase angle, (7 (r — R;, 
Eq. 12.3 and 12.4 can be written as 
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For the present analysis, the aerodynamic force coefficients in Eq. 10 and Eq. 12 are obtained 
using the Unear unsteady aerodynamic model of Ref. 3 and the LNFLX2D of Ref. 19. 

Substituting Eq. 7 and Eq. 1 1 in Eq. 5, equation 5 can be written as 


N 

[M] k {fj} k +[K] k {ij} k = co 2 ^[A] kn {Ti} n + (0 2 {G} k 


n = 1 


(13) 


where 

UE p = b k X< 2 m \trtLi*, If 1 m,p (i4) 


NAS A/TM — 2002-2 1 1499 


7 


It should be noted that Eq. 14 permits the use of mode shapes that differ from blade to blade. 

By writing the above equation for all blades, the equations of motion for the cascade can be 
written as 


[^ ? ]{^} + [^ ? ]{X} = co 2 [A]{X} + q) 2 {G} 


where 


( 15 ) 
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{*}= 



NPx 1 


(16.5) 


The individual blade matrices, [M] k is of size PxP, and consists of the mass matrix of the k* 
blade, with elements given as 

M u = mi, i=l,P (17.1) 

where m, is the i th modal mass. Similarly, [K] k is of size PxP, and consists of the stiffness matrix 
of the k 1 ” blade, with elements given as 

K,, = M a *07(1 + 2i£,), i=l,P ^ 17 - 2 ) 

where to. is the natural frequency of the i' h mode and Q is structural damping ratio of the i mode. 


Each element of [aT and {G}* are given by summing over M strips as 
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Assuming the solution for Eq. 15 is of the form 


{x}=fxy* 


(18.1) 
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(19) 
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and dividing both sides by an assumed frequency,^ , and rearranging, the Eq. 15 can be written 
as: 


[M- rfe]F}=rtc} 


where 


( 20 ) 



(21.1) 




( 21 . 2 ) 


For a stability calculation (flutter), the motion-independent forces {G} are set to zero and the 
eigenvalue problem is obtained in the standard form: 


[M-y[0]]{x}={o} (22) 

The solution of the above eigenvalue problem (22) results in NP complex eigenvalues of the 
form 


( 

. (o 

i — 

V«o 



(23) 


The real part of the eigenvalue (fi) represents the damping ratio, and the imaginary part (v) 
represents the damped frequency; flutter occurs if ju > 0 for any of the eigenvalues. 

The aeroelastic response of the blades induced by wakes is calculated from equation (20) as 

{xHM-xfeirVtG} 


tuned cascade 

For a tuned cascade, in which all the blades have identical structural properties, the interblade 
phase angle modes are uncoupled. The equation of motion can be solved for each interblade 

phase angle, a R for the R* mode, given by Eq. 3. For tuned cascade analysis, Eq, 13 can be 
written as 


[M] k {tj} k +[K] k {ri} k =co 2 [A] R {ri} + (o 2 {G} R (25) 

Since the blades are identical, the same equation is obtained for each blade, and superscript ‘k’ 
can be dropped. The unsteady aerodynamic force coefficients in calculating [A], Eq. 10, are 
given by 
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(26.1) 
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Equation 25 is solved for N different values of the interblade phase angle given by Eq. 3. It is to 
be noted here that the size of the matrices for solution is now reduced to PxP compared to 
(NPxNP) for mistuned case. As before, the equations for the stability (flutter) problem are 
obtained by setting the motion-independent forces to zero. For a given interblade phase angle, 
the solution of the eigenvalue problem results in P complex eigenvalues of the form given by 
Eq. 23 and flutter occurs if p > 0 for any of the eigenvalues. The eigenvalue problem is solved 
for each of the N permissible values. The critical phase angle is identified as the one that results 
in the lowest flutter speed. 

Stability calculation 

The aerodynamic coefficients have to be calculated before the eigenvalue problem can be set up 
and solved. Since the unsteady aerodynamic coefficients depend on the frequency of oscillation, 
it is necessary to assume a frequency co 0 (actual input to the code may be the reduced frequency 
of blade vibration based on chord, k c =ac/U, where U is the free stream velocity ) in advance to 
be able to calculate the aerodynamic coefficients. The aerodynamic coefficients are functions of 
inlet Mach number M, and interblade phase angle <7 r , in addition to cascade geometric 
parameters. In the present study, for a given inlet Mach number, the reduced frequency is varied 
until the real part of one of the eigenvalues p becomes zero while the real parts of the remaining 
eigenvalues are negative. This is repeated for all possible interblade phase angles for a tuned 
cascade. The assumed flutter-reduced frequency k c f and the calculated flutter frequency v f are 
both based on co f . Thus, these two can be combined to eliminate 0) f and the flutter speed is 
obtained, namely, U f = v f c(O 0 /k cf . Since the inlet Mach number is known, this flutter speed 
gives the inlet condition (speed of sound, aoo ) at which the cascade will be neutrally stable for 
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the given Mach number. This procedure can be repeated to obtain a plot of flutter speed versus 
Mach number. Knowing the operating conditions, it is possible to determine whether flutter will 
occur within the operating region and if so, the Mach number and frequency at flutter. It should 
be noted that a mistuned cascade can also be analyzed using this procedure but will require the 
solution of the equations for all phase angles at one time for a given reduced frequency. 

Analysis procedure 

The flutter and forced response analysis in FREPS-LE consists of running five codes, (I) a 
structural dynamic analysis code, such as NASTRAN, ANSYS, etc., (II) 2DSTRIP, (III) 
NPHASE, (IV) LINFLX2D, and (V) 2DASTROP. It is to be noted that ASTROP2 code is a 
combination of 2DSTRIP and 2DASTROP codes. The analysis procedure is explained in five 
steps. In step 1 a vibration analysis is performed for the blade. The output is natural frequencies 
and mode shapes. This output is used by the 2DSTRIP code. In step 2, strips are selected, and 
2DSTRIP is run to calculate relative Mach numbers, sweep angles, stagger angles, chord values, 
and strip widths at these strips. During this run, the three dimensional modal values are also 
interpolated at each strip, and equivalent pitching and plunging modal values are calculated. In 
step 3, a steady aerodynamic solution at these strips is obtained from NPHASE. The steady 
aerodynamic solution is written as a database. Step 4 consists of running LINFLX2D for 
assumed number of reduced frequencies, interblade phase angles, for pitching and plunging 
modes. The analysis is carried out for unit amplitude of vibration for all the strips, and the 
unsteady aerodynamic solutions are stored as a database. In step 5, the 2DASTROP code is 
executed to calculate flutter using the eigenvalue approach and to calculate forced response. The 
number of calculations required for the aeroelastic analysis depend on the number of strips, 
Mach numbers, interblade phase angles, and reduced frequencies as described below. 

Let each blade be divided into NSTRIP number of strips. In general, for a given inlet Mach 
number, the flow conditions, (relative Mach number and angle of attack), and the geometric 
conditions, (gap to chord ratio, stagger angle, and airfoil shape) will be different at these strips. 
Therefore, the number of steady aerodynamic solutions (number of NPHASE runs) required is 

NSTDY = NMACH * NSTRIP, where NMACH is the number of inlet Mach numbers to be 
considered in the study. 

In the case of the unsteady aerodynamic solution, three other parameters, number of modes, 
reduced frequency (k c ) and interblade phase angle (cr) have to be considered. If NMODE is the 
number of modes, NREDF is the number of reduced frequencies, and NSIGMA is the number of 
interblade phase angles, then the number of unsteady solutions for flutter (number of LINFLX2D 
runs) is, NUSTDY = NSTRIP*NMACH * NMODE *NSIGMA * NREDF. It is to be noted that 
NSIGMA is equal to the number of blades of the cascade. The forced response is usually 
calculated for a specific excitation at a fixed interblade phase angle and at a given frequency. 
Therefore, the number of unsteady runs for forced response are reduced to NUSTDY = 
NSTRIP*NMACH, same as for steady solution. 

The main purpose of the present effort is to couple the unsteady aerodynamic solutions from 
LINFLX2D to ASTROP2 code and to validate the procedure. In order to reduce the number of 
calculations required for the aeroelastic analysis, a straight, untapered stator blade at two inlet 
Mach numbers will be considered for the present study. Therefore, for each Mach number, the 
flow and geometric conditions will be same at all strips, reducing the number of steady 
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aerodynamic calculations to NMACH, i.e. NSTDY = 2. For the unsteady run, the number of 
modes is fixed as two, plunging and pitching, i.e. NMODE=2. The number of unsteady runs for 
flutter is given as NUSTDY = 2 *2 *NSIGMA *NREDF = 4*NSIGMA* NREDF. For 
predicting the flutter boundary, calculations have to be carried out for a number of reduced 
frequencies. It is suggested that a LINFLX2D database be prepared for three or four reduced 
frequencies, and then interpolation be used for required frequencies. These reduced frequency 
values can be selected by first running ASTROP2-LE with linear theory. The number of 
unsteady runs for forced response is given by NUSTDY =2. 


Results and Discussion 

Results presented here are meant to demonstrate the state of development of the code and to 
show that the analysis procedure given in the previous sections has been implemented correctly. 
Therefore, calculations are made for a non-rotating cantilevered blade, representing a stator blade 
of a turbomachinery component. Results are presented for an assembly of 12 blades on a rigid 
disk coupled only aerodynamically. Similar geometry is also considered in Ref. 23. 

Both tuned and mistuned cascade analyses are carried out for a flat plate and the tenth standard 
configuration, (designated as CIO, Ref. 19) airfoil cross sections. A 20% alternate mistumng is 
considered for the mistuning analysis. The gap to chord ratio (s/c) is 1.0 and the stagger angle 
(theta) is 45 degrees. Two inlet Mach numbers are considered M=0.7 and 0.8. For the CIO 
airfoil and for M=0.8, the flow is transonic with a shock near the quarter chord. An H-grid of 
141x41 grid is used for the study. There are 80 points on the airfoil, and the inlet boundary is 5 
chords from the airfoil leading edge and the exit boundary is 10 chords from the airfoil trailing 
edge. Steady aerodynamic solutions were obtained for these Mach numbers using NPHASE 
before running LINFLX2D. 

LINFLX2D code: 

To validate the LINFLX2D code, unsteady aerodynamic pressure differences were calculated for 
both airfoil cross sections for a reduced frequency of oscillation based on chord (k c ) of 1.0. The 
unsteady pressure difference is non-dimensionalized by (upper surface pressure- lower surface 
pressure )/(airdensity*U**2*amplitude of oscillation). 

Subsonic Inflow 

Figure 2a shows the steady Mach number distribution obtained for the CIO airfoil cascade for an 
inlet Mach number, M=0.7. The flow to the cascade is at 10 degrees angle of attack. It can be 
seen that for this Mach number the flow is shock free. The unsteady pressure difference 
distribution for an interblade phase angle, <7, of 180 degrees is shown in Fig. 2b. The blades are 
oscillating in pitch about the midchord. Figure 2b shows predictions from Unear theory ( Ret. 4) 
and from the linearized Euler code LINFLX2D for flat plate geometry, from the nonlinear Euler 
(Ref 17) and from LINFLX2D for the CIO airfoil geometry. The unsteady pressure difference 
distribution obtained with LINFLX2D for the flat plate geometry correlates very well with linear 
theory Also, the predictions from LINFLX2D for the CIO geometry correlate well with the 
nonlinear Euler results indicating that the unsteady aerodynamic pressures predicted by 
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LD S IFLX2D are accurate. From this figure the effect of geometry, angle of attack and blade 
thickness can also be seen indicating the importance of including these features in the unsteady 
aerodynamic analysis. 

Transonic Inflow 


For the CIO airfoil cascade geometry, and for an inlet Mach number of 0.8, the flow is transonic 
with a normal shock occurring in each blade passage. The steady Mach number distribution 
obtained from NPHASE is shown in Fig. 3a for a steady angle of attack of 13 degrees. A normal 
shock occurs on the suction surface at about 28% of the chord. Results from the full potential 
solver, Ref. 12, are also included for comparison. Both results agree well, except that the 
nonlinear steady Euler solver, NPHASE, predicts the shock location slightly downstream of that 
predicted by the full potential solver of Ref. 1 2. A shock fitting procedure was used in Ref. 1 2, 
whereas the shock is captured naturally in the present solver. 

The unsteady pressure distribution is shown in Fig. 3b, along with a comparison with linear 
theory (Ref. 4) and the unsteady nonlinear Euler solver. The unsteady results are for pitching 
about mid-chord with k c = 1.0, <7= 180° and amplitude of oscillation, cxo, =2°. As expected, 
linear theory does not show any shock, indicating that unsteady analysis based on linear theory 
will not be accurate for cascade flutter analysis in transonic flow. A fair correlation between 
LINFLX2D and non-linear unsteady Euler results can be seen in Fig. 3b. The differences could 
be due to the fact that for transonic flow the following factors may have more effect than for in 
subsonic flow: (1) grids used or (2) the steady solution on which LINFLX2D based is not a low 
loss solution, which is a requirement for an accurate solution from LINFLX2D. These issues 
need separate study. 


Similar results were obtained for plunging motion for both subsonic and transonic inflows. 

ASTROP2 code 

The ASTROP2 code was validated for flutter prediction in Ref. 9 by calculating the flutter 
boundary of the SR3CX2 advanced propeller. The vibration characteristics of the propeller blade 
were obtained from the finite element structural analysis code NASATRAN. Linear theory of 
Ref. 5 was used for calculating the unsteady aerodynamic forces. For more details on the 
ASTROP2 code, see references 9, 24 and 25. In the present effort, forced response capability 
was added to the ASTROP2 code, and modifications were made for the analysis of mistuned 
cascades. To validate these upgrades, the calculations from ASTROP2 were compared with 
those obtained from MISER code. Ref. 7 in appendix A. The validation was done for both flutter 
and forced response. As can be seen from the results shown in appendix A, the results from the 
current version of ASTROP2 agree well with those obtained from MISER. 

Aeroelastic Calculations 

A 12 blade stator is considered for stability and response calculations. In general, the procedure 
for aeroelastic analysis starts with a database of the stator blade containing the geometry, mode 
shapes, and modal frequencies. Then strips are selected, and 2DSTRIP is run to calculate the 
relative Mach numbers, sweep angles, stagger angles, chord values, and strip widths at these 
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strips. During this run, the three dimensional modal values are interpolated at each strip, and 
equivalent pitching and plunging modal values are obtained. Then the steady and unsteady 
aerodynamic solutions are obtained for these strips, and^used for flutter and response analysis. 

However, the main aim of this paper is to demonstrate the ASTROP2 and LINFLX2D code 
coupling. Therefore, the analysis is carried out for the simplified geometry of straight, untwisted 
stator blades. Mode shapes, frequencies, and modal values are assumed instead of performing a 
detailed analysis. Two assumed modes are used in the analysis. The first mode is a pure torsion 
mode and the second mode is a coupled bending-torsion mode. The natural vibration frequencies 
are 81.376HZ and 148.02 HZ respectively. The blade is divided into 10 strips. Table 1 shows 
the aerodynamic and geometric input parameters at the 10 strips. With identical flow and 
geometric conditions, the only contribution to the generalized aerodynamic force matrix at the 
strip is due to the different modal values at each strip. Table 2 shows the modal values used in 

this study. 

The unsteady force coefficients are calculated using LINFLX2D for harmonic blade vibration in 
plunging and pitching modes for a reduced frequency, k c , of 0.2 based on chord. The reduced 
frequency used to calculate the unsteady aerodynamic coefficients, k c =0.2 corresponds to a 
frequency of 96.07 HZ. The same grids that were used for the unsteady validation and the 
NPHASE steady solution obtained in the previous section were used again. The moment 
coefficients are calculated about the leading edge. The calculated unsteady force coefficients 
along with the modal values at the strips are used by 2DASTROP code to calculate the elements 
of [A] and {G}. The stability is inferred from the eigenvalues. 

Stability calculations 

Figure 4 shows the root locus plot of the eigenvalues for mode 1 calculated for flat plate 
geometry operating at M=0.7. Results obtained from LINFLX2D are compared with those 
obtained using linear theory unsteady aerodynamics. Good correlation is observed. At M=0.7 
the flow field is linear hence, linear theory and linearized Euler are expected to correlate well 
with each other. The linear theory unsteady aerodynamic code was an integral part of 
ASTROP2, and LINFLX2D is the new code that is coupled with ASTROP2. The root locus plot 
correlation shows that the coupling of the LINFLX2D unsteady aerodynamic database to 
ASTROP2 is accurate. It is to be noted that since the first mode is a pure torsional mode, the 
frequency of the aeroelastic system is close to the first natural frequency of the blade. 

To explore the effects of airfoil shape, and transonic flow the analysis was carried out with the 
CIO geometry and for M=0.7 and 0.8. Figure 5a shows the root locus plot for a tuned cascade 
for the first mode at an inflow M of 0.7 for the flat plate and CIO airfoil geometries. For the 
tuned case, it can be seen that the blade is more unstable when the effect of airfoil geometry is 
included. The root locus plot for the first mode for tuned cascade for M=0.8 is shown in Fig. 5b. 
Here the root locus plot shows the effect of airfoil shape, and transonic flow. It is seen that the 
tuned blade is slightly less unstable than in figure 5a when the effects of airfoil geometry, and 
shock are included. 

To include the effect of mistuning, the first mode frequency of the alternate blades is 
increased/decreased by 20% i.e. two adjacent blades have frequency of 1.1 fl and 0.9 fl, where 
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fl is the tuned first mode frequency. The root locus plot for the tuned and mistuned rotors is 
shown in figure 6. The frequencies are plotted on double Y -axis; left Y -axis representing tuned 
cascade frequencies (circles), and right Y-axis representing mistuned cascade frequencies 
(triangles). The open symbols represent flat plate airfoil, and closed symbols represent CIO 
airfoil. 

Figure 6a shows the root locus plot for mistuned cascade for the first mode at an inflow M of 0.7 
for the flat plate and CIO airfoil geometries. For comparison, the root locus plot of the tuned 
cascade is also shown in Fig. 6a. It can be seen that the mistuning increased the spread of the 
blade frequencies for both geometries. It can also be seen that (1) root locus plot has split into 
two parts, one corresponding to high frequency, and the other to low frequency, (2) the rotor with 
the flat plate geometry became stable with the addition of mistuning, (3) the rotor with the CIO 
airfoil geometry also moved towards the stable direction, but the amount of mistuning was not 
sufficient to make the rotor stable. 

The root locus plot for the first mode for M=0.8 for the mistuned cascade is shown in Fig. 6b. 
Here the root locus plot shows the effect of airfoil shape, transonic flow and mistuning. The 
effect of mistuning is same as that for M = 0.7. 

Response Calculations 

The response of the blade to a vortical disturbance is calculated for flat plate geometry and 
compared with that obtained from linear theory. A structural damping ratio of 0.002 is added to 
make the aeroelastic system stable (see Fig. 4-6) and to limit the amplitude at resonance. 
Figure 7 shows the tuned aeroelastic response for R=6 i.e. c =180 degrees for flat plate 
geometry. The unsteady aerodynamic coefficients are obtained at k c = 0.2 for M=0.7. The 
moment coefficients are calculated about the leading edge. The forcing frequency range 
investigated is limited to a small range around the l a mode frequency. For the tuned cascade, the 
response will be entirely in the r = R mode, and all the blades will have equal amplitudes. The 
amplitude of response is a function of the frequency ratio, co / to o . Figure 7 shows the 1 st 
generalized degree of freedom (ql) response obtained using linear theory (Ref. 5) and present 
LINFLX2D code. It can be seen that calculations from linear theory and from LINFLX2D are 
identical indicating that the coupling of LINFLX2D coefficients to ASTROP2 code for response 
calculations is accurate. 

Figure 8 shows the ql response obtained for the CIO airfoil for M=0.7. Comparing the response 
of CIO with that of flat plate, it can be seen that the response has increased about 140%. This is 
due to high steady loading on the airfoil. Figure 8 also shows the response with 20% alternate 
mistuning for both flat plate geometry and CIO geometry. The response has split into two 
resonance peaks, one for odd blades and one for even blades with equal amplitude, and same as 
that without mistuning. This indicates that the effect of aerodynamic damping is very small for 
the example chosen here. 

Figure 9 shows the ql response obtained for the CIO airfoil for M=0.8. Comparing the response 
of CIO with that of flat plate, it can be seen that the response has increased to only about 68%. 
Observing from Fig 8 that the airfoil shape increased the response to about 140%, the decrease in 
the response in Fig 9 can be attributed to the presence of shock. However, the flow being 
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transonic, it is necessary to check the present calculation by running LINFLX2D for a different 
grid. Figure 9 also shows the response calculations with 20% alternate mistuning. Again, the 
response has split into two resonance peaks, on for odd blades and one for even blades with equal 
amplitude, and same as that without mistuning. This again indicates that the effect of 
aerodynamic damping is very small for the example chosen here. 

The response calculations given above were carried out by including the motion dependent 
aerodynamic matrix, [A], which adds the contribution of aerodynamic damping to the response. 
For the case considered here, the calculations were repeated without the contributions of [A]. 
The response showed negligible change in the amplitude. This may be due to low aerodynamic 
damping contribution from [A] (see Figs. 4-6 ) compared to the added structural damping ratio of 
0.002. Note that the calculation of the elements of [A] requires calculation of aerodynamic 
coefficients for all interblade phase angles (equal to the number of blades) where as calculation 
of the elements of {G} requires the aerodynamic forcing coefficient for only one phase angle. 

Computational times 

The unsteady aerodynamic solution times varied from 20 minutes to 75 minutes on an SGI 
workstation. This is directly related to the number of iterations required for convergence. A 
solution for flat plate geometry took more time than for the CIO airfoil geometry. 

Concluding Remarks 


The transonic flow unsteady linearized Euler aerodynamic solver, LINFLX2D, has been 
successfully coupled with the ASTROP2 aeroelastic analysis code. The resulting code, 
ASTROP2-LE, is validated by comparing predictions from linear theory for flat plate geometry. 
Comparison was done for both flutter and forced response. Results were also presented for 
cascades in subsonic and transonic flow for a standard cascade section known as CIO geometry. 
In addition the code is validated by comparing results for MISER code. It is noted that the 
number of LINFLX2D solutions required is directly related to number of strips. Care has to be 
taken in selecting the number of strips to reduce the number of calculations. The following were 
observed during the study: (1) The steady loading due to the airfoil shape and angle of attack 
destabilized the cascade for the Mach numbers considered, (2) The steady loading due to the 
airfoil shape and angle of attack increased the blade response for subsonic flow and decreased 
response for transonic flow compared to that for a flat plate, (3) for the cascade geometry 
considered here, even 20% mistuning has very small effect on flutter and negligible effect on 

response. 
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APPENDIX A 


Comparison of ASTROP2-LE calculations with MISER 


In the present research effort, forced response capability was added to the ASTROP2-LE code, 
and modifications were made for the analysis of mistimed cascades. To validate these upgrades, 
the calculations from ASTROP2-LE were compared in this appendix with those obtained from 
MISER code, Ref. 8. The validation was carried out for both flutter and forced response. 

A uniform blade is considered with a chord of 2.0 units and a span of 1.0 unit. Since the main 
aim is to validate the calculation of the aerodynamic matrices, the mass and stiffness matrices 
were taken from MISER, and used in ASTROP2-LE. Two modes are used in the calculations, 
with first mode as (h,alfa) = (1. 0,0.0) and the second mode as {h, alfa) = (0.0, 1.0) at all strips. 
The Mach number is 0.7. The gap to chord ratio is 1.0 and the stagger angle is 45 degrees. The 
linear theory of Ref. 5 is used to calculate the elements of the aerodynamic matrices. Table A.l 
shows the elements of the mass and stiffness matrices. 

Table A. 2 shows the comparison of the elements of the aerodynamic and force matrices, 
stability, and response predictions for a tuned cascade. A reduced frequency of 0.5 based on 
semichord is (k b ) used. The interblade phase angle is 180.0 degrees. The blades are oscillating 

about the leading edge. The table shows results from ASTROP2-LE and from MISER code. As 
seen from table A.2, the ASTROP2-LE calculations compare very well with those of MISER, 
validating the ASTROP2-LE modifications for flutter and forced response of tuned cascades. 

To validate the code for mistuned cascade analysis, a 12 blade cascade is considered. The blades 
are assumed to have alternate frequency mistuning of 20% in second mode i.e. the frequency of 
every alternate blade is 0.9f 1 and 1 . If 1 where f 1 is tuned blade frequency of the second mode. 

Figure A1 shows the root locus plot, real part of the eigenvalue versus imaginary part of the 
eigenvalue, for 12 bladed cascade with and without mistuning. Second mode is plotted. The 
reduced frequency k b is 0.5, and the blades are pitching about the leading edge. As can be seen 
ASTROP2 predictions compare exactly with MISER calculations. 

Figure A2 shows the amplitude of the second mode plotted with varying frequency ratio. A 
structural damping ratio of 0.002 was used to limit the amplitude of the blades at resonance. The 
interblade phase angle of the forcing function is 180 degrees corresponding to an engine order 
excitation of 6. The reduced frequency k b is 0.5, and the blades are pitching about the leading 
edge. Again, it can be seen that the prediction from ASTROP2 match exactly with MISER 
calculations. 

This validates the ASTROP2 code and the modifications. 
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Table 1: Aerodynamic Input Parameters at the Strips 


Atmospheric pressure (psi)= 1 3. 1023; Speed of sound (fps)-l 130.0 


Strip 

Index 

Stagger 

Angle 

(Degrees) 

Chord Length 
(Inches) 

Gap/Chord 

Ratio 

Radius 

(Inches) 

Strip 

width 

(inches) 

1 

45.0 

3.145 

1.0 

5.407 

0.1330 

2 

45.0 

3.145 

1.0 

5.540 

0.1335 

3 

45.0 

3.145 

1.0 

5.674 

0.1335 

4 

45.0 

3.145 

1.0 

5.807 

0.1335 

5 

45.0 

3.145 

1.0 

5.941 

0.1335 

6 

45.0 

3.145 

1.0 

6.074 

0.1335 

7 

45.0 

3.145 

1.0 

6.208 

0.1335 

8 

45.0 

3.145 

1.0 

6.341 

0.1335 

9 

45.0 

3.145 

1.0 

6.475 

0.1335 

10 

45.0 

3.145 

1.0 

6.608 

0.0669 
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Table 2: Mode shapes and frequencies used in the study 


Table 2.1: Mode 1: natural frequency^ 1.376 HZ 


Strip Index 

Bending 

Torsion 

1 

0.0 

0.328 

2 

0.0 

0.444 

3 

0.0 

0.584 

4 

0.0 

0.627 

5 

0.0 

0.693 

6 

0.0 

0.742 

7 

0.0 

0.793 

8 

0.0 

0.889 

9 

0.0 

0.937 

10 

0.0 

1.000 


Table 2.2: Mode 2: natural frequency- 148.02 HZ 


Strip Index 

Bending 

Torsion 

1 

-0.082 

0.329 

2 

-0.087 

0.446 

3 

-0.070 

0.586 

4 

-0.058 

0.630 

5 

-0.034 

0.696 

6 

0.0 

0.745 

7 

0.076 

0.796 

8 

0.203 

0.892 

9 

0.358 

0.940 

10 

0.502 

1.003 
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Table A. 1 


Mass and Stiffness matrices for a tuned cascade 


Mass matrix 

Stiffness matrix 

258.5 0.0 

0.0 86.1815 

32.9456 0.0 

0.0 86.1815 


Table A.2 


Comparison of MISER and ASTR0P2 calculations for a tuned cascade 


Calculation 

Code 

Aerodynamic matrix, A 

Force matrix, G 

Flutter frequency and 
damping 

Response 

magnitude 

ASTROP2 

11 (-1.7174,-4.7294) 

12 (-12.902,-2.6919) 

21 (0.10711,-3.2800) 

22 (-7.3775,-5.9327) 

11 (-4.2854,7.4543) 
21 (-3.5772,4.2408) 

(0.35816,-0.00325) 

(1.04310,-0.04040) 

0.02553 

0.23416 

MISER 

11 (-1.7174,-4.7294) 

12 (-12.902,-2.6919) 

21 (0.10711,-3.2800) 

22 (-7.3775,-5.9327) 

11 (-4.2854,7.4543) 
21 (-3.5772,4.2408) 

(0.35816,-0.00325) 

(1.04310,-0.04040) 

0.02553 

0.23416 
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X 

Figure la: ASTROP2 coordinate system for a rotating blade 



section A-A 


Figure lb: Section A-A showing rigid pitching (a) and plunging (h) motions for the strip 

(reference axis =leading edge) 
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Figure 2a: Steady Mach number distribution for a subsonic cascade, CIO airfoil, gap to chord 
ratio, s/c=1.0, stagger angle, 9 =45 degrees, inlet Mach number, M =0.7, angle of attack, / =10 
degrees. 



distance along chord, x/c 

Figure 2b: Unsteady pressure difference distribution for a subsonic cascade, pitching about 
midchord, s/c= 1.0, 0=45 deg., M=0.7, k=1.0, interblade phase angle, 0^180 deg., amplitude of 
oscillation, oCq =2 deg., i=0 degrees for flat plate and i=10 degrees for CIO airfoil 
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Figure 3a: Steady Mach number distribution for a transonic cascade; CIO configuration, 
s/c = 1.0, 6=45°,M = 0.8, / = 1 3°. 



Figure 3b: Unsteady pressure difference distribution for a transonic cascade; CIO configuration, 
pitching about mid-chord; parameters as in Fig. 3a, kb = 0.5, < 7 = 180°, Oo = 2°. 
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Figure 4: Root locus plot for 12 blade tuned cascade, first mode, flat plate geometry, kb - 0.1, 
structural damping ratio = 0.0, M=0.7 
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Figure 6a: Root locus plot for al2 blade mistuned cascade, first mode, 20% alternate mistumng; 
circles: tuned cascade; triangles: mistuned cascade; open symbol: flat plate; closed symbol. CIO 
airfoil, kb = 0.1, structural damping ratio =0.0, M=0.7; note: left Y-axis refers to tuned cascade 
frequencies, and right Y-axis refers to mistuned cascade frequencies. 
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Figure 6b: Root locus plot for al2 blade mistuned cascade, first mode, 20% alternate mistuning; 
circles: tuned cascade; triangles: mistuned cascade; open symbol: flat plate; closed symbol: CIO 
airfoil, k}y — 0.1, structural damping ratio =0.0', M=0.8; note: left Y-axis refers to tuned cascade 
frequencies, and right Y-axis refers to mistuned cascade frequencies. 
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Figure 7: First mode response for a tuned cascade, flat plate geometry, structural damping 
=0.002, kb = 0.1, N=12, flat plate geometry, M=0.7 
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Figure 8: First mode response for a tuned cascade, flat plate geometry, structural damping ratio 
=0.002, kb = 0.1, N=12, flat plate and CIO geometries, M=0.7 
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Figure 9: First mode response for a tuned cascade, flat plate geometry, structural damping 
=0.002, kb = 0.1, N=12, flat plate and CIO geometries, M=0.8 
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: Al. root locus plot for the second mode; 12 blades, s/c=1.0, stagger=45. 
Pitching about leading edge, kb=0.5, M=0.7, flat plate 
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Figure A2. forced response, R=6(sigma=180 degrees) 
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